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125?.  Relaxation^metiiods . We  saw  that  the  application  of  the  method 
of  finite  differences  may  lead  to  the  solution  of  certain  systems  of 
linear  algebraic  equations.  Such  systems  frequently  arise  in  applications, 
and  a technique  for  solving  them,  bearing  the  name  of  the  relaxation 
method,  is  the  subject  of  this  section. 

Let  us  begin  by  considering  a very  simple  example,  the  system 
(l25.l),  which  corresponds  to  the  problem  of  static  equilibrium  of  ten 
equal  masses,  equally  spaced  on  a light  string  under  a uniform  tension. 

(See  Fig.  71,  in  which  the  sag  is  exaggerated.)  By  symmetry  we  asstime 
Ug  ■ u^,  u?  “ u^,  etc.  This  example,  together  with  other  material  in 
this  section,  is  adapted  from  G.  E.  Forsythe's  chapter  in  E.  F.  Beckeribach 
(editor).  Modern  Mathematics  for  the  Engineer,  The  McGraw-Hill  Book 
Company,  1955.* 
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the  sponsorship  of  the  Office  of  Naval  Research,  project  number 
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The  desired  equilibrium  position  corresponds  to  the  u^,  •••,  u^ 
obtained  by  solving  (12% .l)  with  r^  * • ••  = r^  **  0.  A direct  approach 
to  solving  the  problem  would  be  to  solve  (12%. l)  by  the  systematic 
elimination  of  unknowns,  in  the  fashion  of  high-school  algebra.  The 
elimination  might  alternatively  be  expressed  in  terms  of  determinants. 
While  these  direct  methods  would  solve  (12%. l)  fairly  readily,  they 
would  be  very  complicated  indeed  for  a larger  system  like  that 
corresponding  to  Figure  72. 

An  alternate  solution  of  (121?.!)  has  proved  very  popular  among 
engineers  and  computers.  We  make  a first  guess  of  the  u*s,  say 
ux  = U,  U2  « 8,  u3  - 12,  uu  * lit,  u^  ■ I1?.  From  (12%. l)  we  find  that 


-1,  r) 


r^  ■ 0.  Since  these  residuals  r^  are  not 


all  zero,  we  will  improve  the  trial  solution  by  reducing  them.  We 
pick  one  of  the  numerically  largest  r^,  say  r^,  and  bring  it  to  zero  by 
an  appropriate  change  Au^  in  u^  alone.  It  is  clear  from  (^^.l)  that 
a unit  increase  in  u^  (Au^  * +l)  would  cause  changes  only  in  r^, 
and  r^,  and  these  changes  would  be  m -1,  Ar^  ■ +2,  Ar^  * -1.  To 

make  r^  ■ 0 calls  for  Ar^  ■ -1,  which  we  bring  about  by  selecting 
Au^  ■ - ,%  . As  by-products  we  have  Arg  ■ Ar^  ■ Accumulating  the 

r's  and  Ar’s,  we  find  the  residuals  r^  = -1.0,  * -0.1?,  r^  ■ 0.0, 

r^  ■ O.v,  r^  - 0.0.  There  is  now  a single  numerically  largest 
residual,  r^,  and  we  proceed  to  "liquidate"  it  by  selecting  Au^  ■ + .E>. 
Next  time  Aug  ■ + .?,  etc.  Eight  steps  of  this  process  are  summarized 
in  Table  1.  An  experienced  computer  goes  very  rapidly,  calculating 
mentally  and  recording  a residual  only  when  it  changes. 
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Table  1.  Relaxation  Solution  of  Equations  (125.1) 
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True  solution 
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Southwell  thinks  of  the  r's  as  negatives  of  constraining  forces 
actually  applied  to  the  weights  to  keep  the  system  in  equilibrium  with 
the  current  displacements.  Each  step  of  the  above  calculation  is  then 
thought  of  as  a relaxation  of  one  of  these  external  constraints.  Hence 
Southwell's  name  for  the  process — relaxation. 

At  the  bottom  of  Table  1 are  cumulated  the  current  values  of  the  u's. 
E.g.,  u^  “ U + .£  + .3  ■ h.8.  Re-calculation  of  the  residuals  then 
confirms  the  computation  so  far.  In  these  eight  mental  steps  the 
maximum  error  of  the  u's  has  been  reduced  from  1.0  to  0.3.  Further 
computing  would  improve  the  u's  at  a comparable  rate,  and  it  would 
not  take  long  to  achieve  ordinary  engineering  accuracy, 

there  are  many  tricks  used  by  relaxers.  One,  illustrated  in  Table 
1,  is  to  work  to  one  significant  digit  only,  and  not  to  complicate  the 
numbers  by  introducing  overprecise  corrections  like  Au^  - -,2S.  Thus 
residuals  are  liquidated  only  in  the  most  significant  digit.  More 
precision  comes  automatically  in  later  steps.  Other  tricks  can  be 
used  to  accelerate  the  convergence  of  the  u's  to  the  correct  answers. 

Such  acceleration  is  nearly  always  essential  to  solving  a problem  of 
any  magnitude. 

A great  time-saver  in  engineering  practice  is  not  to  draw  up 
anything  like  Table  1,  but  instead  to  use  a working  drawing  of  the 

model  as  a computing  sheet.  The  values  of  A^  and  r^  can  be  recorded 
on  the  drawing. 

Relaxation  is  really  fun  for  a computer,  for  several  reasons  I 
(l)  seeing  the  partial  answer  evolve  lends  a purpose  to  each  step,  and 
combats  the  usual  tedium  of  day-long  computing*  (2)  one's  d telligence 


is  continually  challenged  by  the  possibility  of  improving  the  speed  of 
convergence}  (3)  one  need  never  waste  much  time  in  erroneous  computing, 
as  is  possible  in  elimination* 

There  are  many  variations  of  relaxation  methods.  They  all  deal 
with  solving  systems  of  equations,  usually  linear,  and  they  share  these 
essential  properties i (i)  for  any  trial  solution  there  is  a measure  of 
the  error  in  each  of  the  equations}  (ii)  for  each  unsatisfied  equation 
there  is  a separate  formula  for  improving  the  trial  solution}  (iii)  one 
calculates  at  each  step  with  the  equation  whose  error  is  largest. 

The  relaxation  method  was  originally  devised  for  pe ncil-and-paper 
computing,  without  a keyboard  calculator,  and  is  ideally  adapted  to 
such  work.  Tt  is  reasonably  adaptable  to  keyboard  calculators,  but 
here  it  seems  to  lose  some  of  its  relative  superiority  over  other 

methods.  For  automatic  digital  computers,  see  below. 

2 

The  relaxation  method  seems  to  date  from  Gauss,  who  used  and 
recommended  the  basic  method  and  many  of  the  standard  tricks.  Seidel^ 
proved  it  would  converge  for  linear  systems  w.th  positive  definite 
matrices.  In  the  thirties  Southwell^  rediscovered  Gauss’s  method, 
and  named  it.  He  and  his  school  have  developed  the  method  and 
brought  its  wide  applicability  to  the  attention  of  engineers  and 
scientists  everywhere.  The  method  has  proved  especially  suited  to  the 
analysis  of  complicated  redundant,  pin- jointed  frameworks.  These 
have  equations  like  (12? .1),  with  more  involved  coefficients. 
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As  a digression  it  should  be  noted  that  one-dimensional  problems 
like  that  of  Figure  71  should  not  in  practice  be  solved  by  relaxation. 
For  they  are  easily  solved  by  a fast  trial- and- recurs ion  scheme,  even 
for  variable  masses.  For  example,  in  (l)  change  the  last  equation  to 
read 


-1 


* 2u?  - u6  “ r^ 


- 0 


Fix  r^  ■ ••• 


r^  ■ 0.  By  symmetry  we  should  have  Ug  - u^  ■ 0.  Try 

u^  “3.  From  the  modified  equations  (l),  we  find  successively  up^  ■ ?, 

^ (l)  - t’l)  (l)  (l) 

■ 6,  u^  ' ■ ?,  Ug  ' -3,  whence  Ug  ' - u^  ' ■ -2.  Now 


3 uU 

try  up^ . We  find  successively  up^ 


(2) 

u6 


21  u(2)  - u(2) 

Z1»  u6  u? 


11,  u<2>  - 1?,  u[2)  - 18,  u<2>  - 20, 


+1.  Interpolating  linearly  between  up^  and  up> 


to  make  Ug  - u^  zero,  we  get  u^  - whence  the  true  solution  is 
obtained  recursively. 

The  point  of  Table  1 was  to  show  the  technique  of  relaxation  in 
a simple  setting.  The  practical  applications  of  the  method  begin  with 
two-dimensional  problems— like  trusses.  Or,  in  closer  relation  to 
Figure  71,  suppose  one  had  an  L-shaped  network  of  light  strings  with 
21  weights  on  it  in  a horizontal  plane  (Figure  72).  How  might  one 
calculate  the  equilibrium  position  of  the  weights  under  large  tension 
and  under  gravity?  With  the  same  assumptions  as  above,  we  find  21 
equations  like 

(12?  .2)  -1  - Ug  - Ug  + Uuy  - Ug  - ul2  “ r^  ■ 0 . 

It  would  be  a tedious  computation  at  a desk  to  solve  these  by  elimination, 
and  no  simple  recursive  scheme  works  here.  A practical  answer  is 


«£> 
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relaxation,  which  works  numerically  about  like  Table  1,  although  with 
slower  convergence.  The  same  tricks  apply  as  before. 

An  essential  feature  of  systems  of  equations  like  (12?.2)  is  that 
most  of  the  coefficients  of  the  u’s  are  zero.  It  is  this  feature  of 
such  systems  which  makes  relaxation  a possible  pencil-and-paper  method 
of  solving  them. 

The  matrix  A of  coefficients  a^  of  the  system  (12?. 1)  or  of  the 
system  of  all  21  equations  like  (12?.2)  has  two  properties  which  will 
prove  very  important  in  our  further  discussion.  First, 

(12? .3)  A is  symmetric  and  positive  definite. 

The  second  property;  (12?  .U),  concerns  the  geometry  of  the  connecting 
strings  in  Figures  71  and  72.  Note  that  the  weights  in  the  figures  have 
been  drawn  in  two  colors*  black  and  white.  Note  that  each  string  connects 
weights  of  opposite  color.  Hence,  in  the  equations  the  subscripts  i 
of  the  unknowns  u^  can  be  divided  into  two  groups  B,  W (by  color),  so 
that 

(12? .U)  atJ  = 0,  for  i in  B,  j in  B (i  / j)  and  for  i in  W,  J in  W (i  / j) . 

Another  way  of  expressing  (12?. U)  comes  from  reordering  the 
unknowns  ui  and  the  corresponding  equations  so  that  the  "blacks" 
entirely  precede  the  "whites" . Then  the  matrix  takes  the  schematic 
form  of  Figure  73,  where  the  circles  denote  zeros,  the  small  crosses 
denote  nonzero  numbers,  and  the  large  crosses  denote  submatrices  of 
zero  and  nonzero  elements.  Any  system  of  linear  equations  satisfying 
(12?. U)  is  said  by  Young"*  to  have  Property  (A). 
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We  note  in  passing  that  the  first  boundary-value  problem  for  any 

second-order  self-adjoint  partial  differential  equation  lacking  a term 
2 

in  o u/}*)y  leads  to  a symmetric  linear  system  with  Property  (A),  when 
difference  equations  are  suitably  introduced.  If  the  partial  differential 
equation  is  elliptic,  then  (125.3)  holds. 

Before  we  discuss  methods  suitable  for  electronic  computers  it  will 
be  convenient  to  introduce  another  method  for  solving  a linear  system. 

The  most  general  system  of  n linear  algebraic  equations  in  n unknowns 
can  be  written  in  the  form 


(125.5) 


o 


(i  = 1,  •••  , n)  . 


For  the  moment  we  do  not  assume  that  the  matrix  satisfies  (125.3)  or 
(125 .U),  but  it  is  essential  that  no  a^  * 0.  Iterative  methods  for 
solving  (125.5)  have  been  popular  since  Gauss's  time,  if  not  longer. 

One  process,  called  the  Seidel  or  Gauss-Seldel  method,  is  the  following. 
One  solves  the  first  equation  (125.5)  for  u^,  using  the  current  values 
of  Ug,  •••  , un#  Then  the  second  equation  is  solved  for  Ug,  using  the 
latest  known  values  of  u^,  u^,  u^,  •••  , un»  And  so  on.  All  the 
equations  (125.5)  are  solved  in  cyclic  order  for  u^,  •••  , un,  always 
with  the  latest  values  of  the  other  unknowns.  In  other  words,  suppose 
ujk^ , •••  , u^  are  known.  One  gets  u|k+^^ , •••  u^^by  successively 
solving  these  n equations « 
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+1)  + Z UJk)  +bi  = 0 (i  “ 2,  , n 

j-T+1 

(k) 

One  hopes  that  the  will  converge  as  k cd  to  the  which  solve 

the  system  (121>.£). 

The  reader  will  note  that  the  Seidel  process  is  closely  related  to 
the  relaxation  process  described  earlier  in  this  section.  The  difference 
is  the  order  in  which  the  equations  (12S.!>)  are  solved.  In  the 
relaxation  method  we  do  the  equations  in  an  order  determined  by  the  size 
°T  the  r^.  In  the  Seidel  process  the  order  is  fixed  and  cyclic. 

Let  us  analyze  the  behavior  of  the  Seidel  process.  In  our 
applications  the  following  theorem,  apparently  first  completely  proved 
by  Schmeidler^  in  19U9,  is  essential* 

the  A symmetric  and  positive  definite,  then  in  the 

Seidel  process  the  u|^ converge  as  k — > oo  to  limits  u^ 

(i  • , n)  solving  the  system  (12S.!?). 

What  happens  when  A is  not  positive  definite?  Or,  in  any  case, 
how  does  the  vector  u^^  approach  u?  If  the  convergence  is  slow,  how 
can  it  be  speeded  up?  The  answers  to  these  three  questions  can  be 
obtained  through  some  use  of  matrix  theory,  as  follows. 

Let  the  matrix  of  coefficients  A in  (l 2?.5>),  assumed  non— singular, 
be  written  as  the  sum  of  three  matrices:  A - D + E + F.  Here  D has  the 
main  diagonal  of  A,  but  is  zero  elsewhere}  E has  the  below— diagonal 
elements  of  A,  but  is  zero  elsewhere}  and  F has  the  above— diagonal 
elements  of  A,  but  is  zero  elsewhere.  Thus,  schematically. 


(12?.6)  ^ 
j-1 


Equations  (12?. 6)  can  be  written  in  the  following  matrix-vector  form: 
(12?. 7)  (D  + E)u(k+l5  + F u<k>  ♦ b - 0 . 

If  u denotes  the  unique  vector  solving  (12?.?),  we  have 
(12?. 8)  (D  + E)  u ♦ F u +'  b ■ 0 . 

Subtracting  (l2?.8)  from  (12?. 7),  and  letting  e^  - u^  - u denote  the 

(k) 

error  of  u , we  have 

(12?, 9)  (D  ♦ E)  e(k+l)  + F e(k)  - 0 . 

Now  since  no  a^  • 0,  the  matrix  D + E has  an  inverse  (D  + E)”\ 

Letting  H denote  the  matrix  -(D  + E)-1F,  we  find  from  (l2?.9)  that 

(12?. 10)  Q(k+^)  . H e(k)  f 

whence 

(12? .11)  e(k)  - Hk  e(0)  . 

Equation  (12?,10)  shows  the  linear  behavior  of  the  Seidel 
iteration  process.  A representation  of  the  error  e^  in  terms  of  the 
initial  error  is  given  by  (12?. 11),  on  which  a complete  error  analysis 
can  be  based.  The  Gauss-Southwell  relaxation  process  is  theoretically 
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more  complicated  just  because  it  has  no  simple  analog  to  (12? .11). 

From  the  theory  of  linear  transformations  we  know  when  and  how 

(k) 

e goes  to  zero.  The  two  determinants!  equations 

(12*?. 12)  |H  - yUl  | - 0 and  | (D  + E^q  + F|  - 0 

have  the  same  n real  or  complex  roots  •••  , If  all  |yu^  | <1, 

then  in  the  Seidel  process  u^  — > u.  If  any  j * 1,  the  Seidel 
process  diverges.  In  principle  this  settles  the  question  of  convergence. 

For  most  matrices  A,  to  each  of  the  roots ju^  of  (12?.12)  there 
corresponds  a unique  vector  y^  such  that  Hy^  - p y^  . That  is, 
the  transformation  H leaves  the  vector  y^^  unchanged  in  direction, 
but  stretches  it  > l)  or  shrinks  it  ( ^ | < 1)  to  the  fraction 

of  its  previous  length.  All  these  vectors  y^  form  an  oblique 
coordinate  system,  in  terms  of  which  we  can  resolve  the  initial  error 
vector  e^t 

e(0,-£ciy(i)  . 

i-1  1 

Assume  |^|  > \y^  J * ...  i lyMn|.  After  repeated  multiplications  by 
H,  the  resulting  vector  Hk  e^  is  approximately  moved  into  the 
direction  corresponding  to  the  root  of  largest  absolute  value  • 

Hence  we  find  that 

(12?.  13)  e(k)  - Hk  e(0)  - Cl/^k  y(l)  , 

and  we  know  how  fast  e^  ->0.  If  Jk  | < 1,  ultimately  each  step 
reduces  the  length  of  e^  to  the  fraction  J^l  of  itself. 


Hence 


If  J^jJ  < 1,  0 along  one  direction,  that  of  y^\ 

u(k^  _*u  along  the  direction  of  y^.  Cases  where  more  than  one 
|^i  | dondna+e  are  more  complicated,  but  can  be  treated  with  similar 
tools. 

Knowing  the  geometric  character  of  the  convergence,  it  is  not 
difficult  to  design  acceleration  processes  to  speed  up  the  convergence 
of  u^  to  u. 

Xs  an  example  of  the  Seidel  process  and  its  convergence,  we  use 
it  to  solve  the  system  (125.1)  with  the  same  start  as  in  Table  1.  We 
have  the  following  iteration* 


v. 


Several  rounds  of  this  are  shown  in  Table  2.  The  residuals  are  not 
shown,  and  one  line  of  the  table  amounts  to  a full  cycle  of  the 
above  algorithm.. 


Table  2,  Seidel  Solution  of  Equations  (I2?.l) 


worst 


In  the  first  column  of  Table  2 is  given  the  ratio  of  the  worst 

error  of  the  u’s  in  the  preceding  and  following  rows  of  the  table. 

By  (12?,13)  this  ratio  converges  to  (A  solution  of  either 

equation  (12? .12)  gives  ^ ■ ,90U?,  ja^  «-  ,3U£5, ^ " °») 

( k) 

approach  of  u to  u is  one-sided  and  very  regular.  It  will  take 
about  22  cycles  to  gain  one  decimal  point  in  accuracy.  Making 
educated  guesses  at  u in  such  problems  is  easy  in  desk  work,  if 
one  knows  (12? ,13)  and  its  analog  when  |jd^J  « | • 


Ihe  present  availability  of  electronic  digital  computing  machines 
makes  it  possible  to  solve  much  larger  problems  than  have  been 
previously  feasible.  Such  machines  carry  out  arithmetic  operations 
at  an  effective  speed  on  the  order  of  10^  times  faster  than  a human 
being  with  a desk  calculator.  Something  like  10^  numbers  of  desk 
calculator  precision  can  be  held  in  a fast— access  "memory"  and  made 
available  as  rapidly  as  the  arithmetic  organ  can  operate.  Something 
like  10^  more  numbers  can  be  held  in  an  intermediate  storage  and 
transferred  to  the  hi^i-speed  memory  in  a few  ndlliseoonds.  Moreover, 
current  developments  will  probably  have  made  the  figures  given  here 
obsolete  before  this  book  is  published. 

Because  of  the  speed  and  capacity  of  such  computers,  many  persons 
want  to  solve  their  problems  on  them.  It  is  pertinent  to  ask.  What 
methods  will  prove  most  feasible  for  the  computers?  While  definitive 
answers  must  await  investigations  as  yet  undone,  certain  indications 
are  now  possible. 

A first  observation  is  that  for  large  problems  of  the  types  of 
(12?. l)  or  (l2?.2),  iterative  methods  are  relatively  attractive, 
for  much  the  same  reasons  as  for  pencil-and-paper  calculation. 

But  the  relaxation  method  as  outlined  in  connection  with  (12?. l) 
has  one  considerable  disadvantage.  The  scanning  of  all  the 
residuals  r^  in  a search  for  max^  |r^  | is  comparatively  time- 
consuming.  In  fact,  while  computing  r^  it  would  take  almost  no 
extra  time  to  solve  the  i-th  equation  for  u^.  But  if  one  solves 
the  i-th  equation  for  u^  (i  ■ 1,  •••  , n),  one  is  actually  carrying 
out  the  Seidel  process,  which  is  accordingly  preferred  in  machine 
calculation  to  conventional  relaxation. 
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A second  observation  is  that  solving  a large  system  (I 25.5)  by 
the  Seidel  method  is  likely  to  be  slow*  To  speed  up  the  solution; 
acceleration  methods  are  needed,  as  indicated  above.  But  accelerations 
involve  new  routines,  new  coding,  and  the  mundane  but  important 
problems  of  storing  or  reading  in  new  codes.  It  is  important  with 
machines  to  reduce  coding  and  operating  to  the  utmost  in  simplicity. 

It  is  the  remarkable  discovery  of  Young  that  for  certain 
problems  a modification  of  the  Seidel  process  will  vastly  speed  up 
the  convergence,  while  scarcely  complicating  the  coding  at  all. 

Where  applicable,  it  thus  eliminates  the  necessity  for  special 
acceleration  routines. 

It  has  long  been  observed  by  relaxers  that  the  Gauss-Southwell 

process  usually  goes  faster  if  one  "overrelaxes"  a little  at  each 

step.  Young  was  therefore  convinced  of  the  value  of  analyzing 

overrelaxation  carefully  in  conjunction  with  the  Seidel  process.. 

He  confines  himself  to  positive  definite  symmetric  matrices  with 

Property  (A)j  see  (125. U) • In  each  step  of  the  Seidel  process 

(125.6)  or  (125.7),  Young  suggests  that  one  first  compute  the 

(k) 

Seidel  value  — call  it  v)  ' — and  then  compute 


4k*1)  - u<k)  ♦ /J(T<k>  - „<k))  . 


This  amounts  to  an  overrelaxation  of  100(  - l)  per  cent.  Young 

asks  which  choice  of  ft  (l  < fl  < 2)  is  best. 

The  analysis  proceeds  much  as  in  the  Seidel  process,  since  this 
systematic  overrelaxation  process  is  also  a linear  one . We  have 


& 
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(12?.1U)  jjjjCk+i)  «.  ^(k)  + ^(k)  + b . o 

and 

■ u(k>  * /3[T(«  - u'«]  - (i  -/! ) U<»  * /?„<k>  . 

We  now  eliminate  from  (12? ,1U).  Since 

j^Of+D  . (].  _ ^ ) Dn^k)  * /5Dv^  , 

we  hare 

(125 .15)  E0(k"1>  ♦ ^ Du^l)  * (l  - 4 j Du<k>  - Pa<k) 

Equation  (12?. 1?)  describes  systematic  overrelaxation.  Just  as 
equation  (12?. 7)  describes  the  Seidel  process.  (Note  that  (12?.1?) 
reduces  to  (12?.7)  for  1.)  The  speed  of  the  convergence  of 

Young's  process  is  measured  by  the  largest  in  modulus  of  the  roots 
^1®  following  deterndnantal  equation,  analogous  to  the  second 
part  of  (12?.12)« 

(12? .16)  |<TE+2D+^l-^ljD+F|-0  . 

For  matrices  satisfying  (12?.3)  and  (12?.U)  and  for  a certain 
ordering  of  the  equations,  it  can  be  shown  that  the  maximum  of  the 
1*1 1 is  least  when  we  choose  /?  - 2(1  + vl  - f1,  where  is  the 
largest  root  of  (12? .12).  Hence  this  defines  the  optimal  amount 
to  overrelax.  Moreover,  for  this  all  |c^  | are  equal. 
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To  illustrate  the  method,  we  show  in  Table  3 the  result  of 
solving  (125.1)  by  systematic  overrelaxation.  Corresponding  to 
» .901*5,  we  take  - 1.528. 

Table  3*  Solution  of  (125.1)  by  Systematic  Overrelaxation 


RatiA 
of  \ 
worst  ) 
errors 

u 

\ 

5 

9 

12 

lli 

15 

u(0> 

U 

8 

12 

Hi 

15 

-3U 

U.76 

9.3U 

12.26 

1U.20 

15.31 

.51 

u<2> 

5.39 

9.32 

12.26 

Hi  .33 

15  .3li 

.51* 

u(3) 

5.039 

9.060 

12.161 

lli  .209 

15.H1O 

.76 

u‘W 

5.025 

9.110 

12.159 

Hi  .119 

15.108 

.7U 

5.071 

9.118 

12 .097 

lli  .093 

15.085 

.52 

5.053 

9.052 

12.059 

ll*  .061 

15.01,8 

.59 

k(7) 

5.012 

9.027 

12 .036 

lli.032 

15.021, 

.67 

„<«> 

5. oil* 

9.021* 

12.021* 

Hi  .020 

15.018 

.58 

u(9> 

5.011 

9. oil* 

12.013 

Hi. 013 

15.010 

.56 

„(10) 

5.001*9 

9.0063 

12.0079 

11*  .0068 

15.0051 
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The  same  ratio  of  worst  errors  is  given  as  in  Table  2.  It  is 
irregular  in  Table  3,  but  will  converge  to  ft  - 1 - .528.  With 
this  value,  one  will  require  only  3.6  cycles  per  decimal  point* 

Thus  the  value  ft  ■ 1.528  is  much  superior  to  the  value  ft  m 1 
of  the  Seidel  method. 

O 

G.  E.  Forsythe  has  coded  for  SWAC°  the  Young  method  for  a simple 
difference  equation, 

(1207) 

u(x  ♦ h,  y)  + u(x,  y + h)  + u(x  - h,  y)  + u(x,  y - h)  - U n(x,  y)  - 0 , 

corresponding  to  the  Laplace  differential  equation.  The  code  will 
accommodate  a network  as  large  as  32  by  128  points.  The  shape  of 
the  boundary  is  inn  ate  rial,  except  that  the  boundary  points  must 
lie  on  the  nodes  of  the  network.  The  following  exaiqple  may 
illustrate  the  usefulness  of  the  method. 

For  a rectangle  of  30  X 68  • 20lj0  interior  unknown  points,  each 
cycle  of  relaxation  takes  8.5  seconds.  For  the  Seidel  process  the 
dominant  eigenvalue^  - 0.99ljl6.  To  reduce  the  error  to  10"6  times 
its  initial  value  by  the  Seidel  process  ( ft  m l)  would  require 
about  2300  cycles,  because  approximately  one  has  ( .99U16)2-300  - 1(T^. 
This  would  take  over  five  hours  on  SWAC.  If,  however,  ft  is  taken 
at  its  optimal  value  of  2(1  + (l  - ,99la6)^)_1  - 1.85802,  then  the 
dominant  eigenvalue  is  <5^  - .85802.  For  this  <5^,  it  requires  about 
90  iterations,  accomplished  in  only  13  minutes,  to  reduce  the  error 
from  1 to  10"^. 
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Some  practice  with  SWAC  convinces  people  it  is  not 

difficult  to  estimate  ft  veil  enough  in  a few  minutes . 

Actual  running  time  to  reduce  the  error  by  a factor  of  10“^  is  on 

the  order  of  20  to  30  minutes,  including  the  time  necessary  to 

determine  ft  approximately.  Similar  experience  is  reported  by 
o 

Young  and  Lerch. 

On  SWAC  the  progress  of  the  calculation  can  be  monitored  by 

observing  the  value  of  lri^  I as  bb®  calculation 

proceeds.  When  ft  ■ 1,  E^)  decreases  mono  tonic  ally  and  smoothly. 

When  ft  is  around  the  optimal  value,  has  comparatively  wild 

fluctuations,  probably  because  of  "beats"  (something  of  this  is 

seen  in  Table  3)  .-  between  the  various  complex  frequencies  of 

equal  magnitude.  It  is  perhaps  also  because  the  operator  carrying 
(y)  (k+ll 

E into  Ev  ' has  a nonlinear  elementary  divisor* 

For  recent  research  on  Young's  and  similar  processes,  see 
Riley^-  and  Sheldon.^  Ihere  is  evidence^  that  systematic 
overrelaxation  is  also  useful  for  some  matrices  not  having 
Property  (A).  The  Seidel  method.  Young's  modification,  and  the 
SWAC  codes  are  adaptable  to  obtaining  the  fundamental  eigenvalue 
of  a matrix  with  Property  (A). 


March  30,  19?? 


20 


FOOTNOTES  FOR  SECTION  12? 

■'"Southwell,  R.  7.i  "Relaxation  Methods  in  Engineering  Science," 
Oxford  University  Press,  191*0,  2?2  pp. 

2 

Gauss,  C.  F.f  Brief  an  Gerling,  26  December  1823,  "Werke,"  Vol.  9, 
pp.  278-281.  [Translated  by  G.  E.  Forsythe,  Math  . Tables  and  Other  Aids 
to  Computation.  Vol.  ? (l9?l),  pp.  2??-2?8.] 

3 

Seidel,  Ludwigt  Debar  ein  Verfahren,  die  Gleiohungen,  auf  velche 
die  Methode  der  kleinsten  Quadrate  fbhrt,  sowie  linelre  Gleichungen 
ftberhaupt , durch  successive  Ahnlherung  aufzulBsen,  Abh.  Math.  Phys. 

IQ.,  Bayrische  Akad. Wlss..  Hftnchen.  Vol.  11  (I87U),  pp.  81-108. 

^Southwell,  R.  V.*  op.  cit.  and  "Relaxation  Methods  in  Theoretical 
Physics,"  Clarendon  Press,  Oxford,  19U6,  21*8  pp. 

Young,  Davidr  Iterative  Methods  for  Solving  Partial  Difference 
Equations  of  Elliptic  Type,  Trans.  Amsr.  Math.  Soc . , Vol.  76  (19?U), 
pp.  92—111.  (Condensation  of  his  19^0  Harvard  thesis.) 

^Schmeidler,  Werner:  "Vortrlge  hber  Detenninanten  und  Matrizen  ndt 
Anwendungen  in  Physik  und  Technik,"  AkademieVerlag,  Berlin,  19U9,  1??  pp. 

Professor  A.  Ostrovski  has  traced  incomplete  proofs  back  to  P.  Pizzetti, 
Atti  della  Reale  Accadenda  dei  Lincei.  Rendiconti  (U)  Vol.  30(l887), 
pp.  230-23?,  288-293. 

7 

Op.  cit.,  footnote  ?.  The  same  suggestion  was  put  forward  a 
little  earlier  in  a special  case  by  Frankel,  Stanley  P.i  Convergence 
Rates  of  Iterative  Treatments  of  Partial  Differential  Equations, 

Math.  Tables  and  Other  Aids  to  Computation.  Vol.  U (19?0),  pp.  6?-7? . 


21 


National  Bureau  of  Standards  Western  Automatic  Computer,  located 
In  the  U.C.L.A.  Department  of  Mathematics. 

9 

Young,  David  and  Lerch,  Francis*  The  Numerical  Solution  of 
Laplace's  "Equation  on  ORDVAC,  BaUisticss  Research  Laboratories  Mem  or  andu  m 
Report  No . 708,  Aberdeen  Proving  Ground,  July  1953,  25  pp. 

10Riley,  J.  D.*  Iteration  procedures  for  the  Dirichlet  difference 
problem.  Math.  Tables  and  Other  Aids  to  Computation.  Vol.  8 (195U), 
pp.  125-131. 


11Sheldon,  J.  W.t  On  the  numerical  solution  of  elliptic  difference 
equations.  Math.  Tables  and  Other  Aids  to  Computation.  Vol.  9 (1955), 


pp.  000-000. 

12 

Charney,  J.  G.  and  Phillips,  N.  A.*  Numerical  integration  of 
the  quasi-geos trophic  equations  for  barotropic  and  simple  baroclinic 
flows,  J.  Meteorology.  Vol.  10  (1953),  pp.  71-99. 


